%% Import data from text file.
% Script for importing data from the following text file:
%
%    /Users/Grace/Desktop/GOM Aug 2010 recalculated/GUST matlab
%    software/S4882 8FH4 core B_results.csv
%
% To extend the code to different selected data or a different text file,
% generate a function instead of a script.

% Auto-generated by MATLAB on 2014/05/20 12:48:47

%% Initialize variables.
filename = [
    'S4873 8N04 core B_results.csv';
    'S4874 8NO2 core B_results.csv';  %note:  Removed "D" from Station names 
    'S4876 8H02 core B_results.csv';
    'S4877 8H04 core B_results.csv';
    'S4878 8BH2 core B_results.csv';
    'S4879 8BH4 core B_results.csv';
    'S4881 8FH3 core B_results.csv';
    'S4882 8FH4 core B_results.csv']

station= [
        'H7-1';
        'H7-2';
        'E4-1';
        'E4-2';
        'D5-1';
        'D5-2';
        'A6-1';
        'A6-2']

    [r,c]=size(filename)
for i=1:r
delimiter = ',';
startRow = 2;

%% Format string for each line of text:
%   column1: double (%f)
%	column2: double (%f)
%   column3: double (%f)
%	column4: double (%f)
%   column5: double (%f)
%	column6: double (%f)
%   column7: double (%f)
%	column8: double (%f)
%   column9: double (%f)
%	column10: double (%f)
% For more information, see the TEXTSCAN documentation.
formatSpec = '%f%f%f%f%f%f%f%f%f%f%[^\n\r]';

%% Open the text file.
fileID = fopen(filename(i,:),'r');

%% Read columns of data according to format string.
% This call is based on the structure of the file used to generate this
% code. If an error occurs for a different file, try regenerating the code
% from the Import Tool.
dataArray = textscan(fileID, formatSpec, 'Delimiter', delimiter, 'EmptyValue' ,NaN,'HeaderLines' ,startRow-1, 'ReturnOnError', false);

%% Close the text file.
fclose(fileID);

%% Post processing for unimportable data.
% No unimportable data rules were applied during the import, so no post
% processing code is included. To generate code which works for
% unimportable data, select unimportable cells in a file and regenerate the
% script.

%% Allocate imported array to column variable names
step1(:,i) = dataArray{:, 1};
tau_applied(:,i) = dataArray{:, 2};
time(:,i) = dataArray{:, 3};
e_cb(:,i) = dataArray{:, 4};
e_ce(:,i) = dataArray{:, 5};
lambda(:,i) = dataArray{:, 6};
tau_ce(:,i) = dataArray{:, 7};
tau_e (:,i)= dataArray{:, 8};
E_0 (:,i)= dataArray{:, 9};
M(:,i) = dataArray{:, 10};

%% Clear temporary variables
clearvars delimiter startRow formatSpec fileID dataArray ans;
end

% Create figure
figure1 = figure;

% Create axes
axes1 = axes('Parent',figure1,'YDir','reverse');
box(axes1,'on');
hold(axes1,'all');

plot(tau_applied(1:6,1),e_ce(1:6,1),'b-*','markersize', 8)
plot(tau_applied(1:6,2),e_ce(1:6,2),'b--*','markersize', 8)
plot(tau_applied(1:6,3),e_ce(1:6,3),'r-p','markersize', 8)
plot(tau_applied(1:6,4),e_ce(1:6,4),'r--p','markersize', 8)
plot(tau_applied(1:6,5),e_ce(1:6,5),'g-d','markersize', 8)
plot(tau_applied(1:6,6),e_ce(1:6,6),'g--d','markersize', 8)
plot(tau_applied(1:6,7),e_ce(1:6,7),'k-o','markersize', 8)
plot(tau_applied(1:6,8),e_ce(1:6,8),'k--o','markersize', 8)



legend('H7-1','H7-2','E4-1','E4-2','D5-1','D5-2','A6-1','A6-2')
xlabel('Tau Applied  (Pa)')
ylabel('Eroded Mass  (Kg/m{2})')

print (figure(1), '-dpdf', 'All_tau_applied_VS_eroded_M')

% Create figure 2
figure2 = figure;

% Create axes
axes2 = axes('Parent',figure2,'YDir','reverse');
box(axes2,'on');
hold(axes2,'all');



plot(tau_ce(1:6,1),e_ce(1:6,1),'b-*','markersize', 8)
plot(tau_ce(1:6,2),e_ce(1:6,2),'b--*','markersize', 8)
plot(tau_ce(1:6,3),e_ce(1:6,3),'r-p','markersize', 8)
plot(tau_ce(1:6,4),e_ce(1:6,4),'r--p','markersize', 8)
plot(tau_ce(1:6,5),e_ce(1:6,5),'g-d','markersize', 8)
plot(tau_ce(1:6,6),e_ce(1:6,6),'g--d','markersize', 8)
plot(tau_ce(1:6,7),e_ce(1:6,7),'k-o','markersize', 8)
plot(tau_ce(1:6,8),e_ce(1:6,8),'k--o','markersize', 8)


legend('H7-1','H7-2','E4-1','D5-1','D5-2','A6-1','A6-2')
xlabel('Tau_{cr}  (Pa)')
ylabel('Eroded Mass  (Kg/m^{2})')

print (figure(2), '-dpdf', 'All_tau_crit_VS_eroded_M')

%regression curve tau_ce vs eroded mass for 2014 paper with Kevin Briggs

H7=[tau_ce(3:6,1),e_ce(3:6,1)]
H7(5:8,:)=[tau_ce(3:6,2),e_ce(3:6,2)]

E4=[tau_ce(3:6,3),e_ce(3:6,3)]
E4(5:8,:)=[tau_ce(3:6,4),e_ce(3:6,4)]

D5=[tau_ce(3:6,5),e_ce(3:6,5)]
D5(5:8,:)=[tau_ce(3:6,6),e_ce(3:6,6)]

A6=[tau_ce(3:6,7),e_ce(3:6,7)]
A6(5:8,:)=[tau_ce(3:6,8),e_ce(3:6,8)]

save ('tau_vs_er_mass_H7.asc' H7 '-ascii')


% Create figure
figure3 = figure;

% Create axes
axes1 = axes('Parent',figure3);
box(axes1,'on');
hold(axes1,'all');

plot(A6(:,1),A6(:,2),'ko','markersize',8)
hold on
plot(D5(:,1),D5(:,2),'gd','markersize',8)
plot(E4(:,1),E4(:,2),'rp','markersize',8)
plot(H7(:,1),H7(:,2),'*','markersize',8)


legend('A6','D5','E4','H7')
xlabel('Tau_{cr}  (Pa)', 'fontsize', 16)
ylabel('Eroded Mass  (Kg/m^{2})', 'fontsize', 16, 'rotation', 90)

% Find x values for plotting the fit based on xlim
axesLimits1 = xlim(axes1);
xplot1 = linspace(axesLimits1(1), axesLimits1(2));


fitResultsH7 = polyfit(H7(:,1),H7(:,2), 2);
% Evaluate polynomial
yplot1 = polyval(fitResultsH7, xplot1);
% Plot the fit
fitLine1 = plot(xplot1,yplot1,'--','linewidth', 1.5, 'DisplayName','   quadratic','Parent',axes1,...
    'Tag','quadratic',...
    'Color',[0 0 1]);


fitResultsE4 = polyfit(E4(:,1),E4(:,2), 2);
% Evaluate polynomial
yplot2 = polyval(fitResultsE4, xplot1);
% Plot the fit
fitLine2 = plot(xplot1,yplot2,':','linewidth', 1.5, 'DisplayName','   quadratic','Parent',axes1,...
    'Tag','quadratic',...
    'Color',[1 0 0]);

fitResultsD5 = polyfit(D5(:,1),D5(:,2), 2);
% Evaluate polynomial
yplot3 = polyval(fitResultsD5, xplot1);
% Plot the fit
fitLine3 = plot(xplot1,yplot3,'-.','linewidth', 1.5,'DisplayName','   quadratic','Parent',axes1,...
    'Tag','quadratic',...
    'Color',[0 1 0]);

fitResultsA6 = polyfit(A6(:,1),A6(:,2), 2);
% Evaluate polynomial
yplot4 = polyval(fitResultsA6, xplot1);
% Plot the fit
fitLine4 = plot(xplot1,yplot4,'linewidth', 1.5,'DisplayName','   quadratic','Parent',axes1,...
    'Tag','quadratic',...
    'Color',[0 0 0]);

%create a better figure for print

set(gca, 'Box', 'off' ); % here gca means get current axis
 set(gca, 'TickDir', 'out', 'XTick', [0:0.1:0.5], 'YTick', [0:.05:0.2]);  %set the tick marks so they are outside the box


print (figure(3), '-dpdf', 'All_tau_crit_VS_eroded_M_ASregression_updated')
print (figure(1), '-deps', 'All_tau_crit_VS_eroded_M_ASregression_updated')
print (figure(1), '-djpeg', 'All_tau_crit_VS_eroded_M_ASregression_updated')
